import numpy as np 
import numpy.random as npr
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
from scipy.optimize import leastsq,least_squares

def logistic4(x,A,B,C,D):
    '''4PL logistic equation'''
    return ((A-D)/(1.0+((x/C)**B))) + D

def logistic5(x,p1,p2,p3,p4,p5):
    c = 2*p3*p5/abs(p3+p5)
    f = 1/(1 + np.exp(-c*(p4-x)))
    g = np.exp(p3*(p4-x))
    h = np.exp(p5*(p4-x))
    return p1 + (p2/(1 + f*g + (1-f)*h))

def residuals(p,y,x):
    A,B,C,D = p
    err = y - logistic4(x,A,B,C,D)
    return err

def residuals5(p,y,x):
    p1,p2,p3,p4,p5 = p
    err = y - logistic5(x,p1,p2,p3,p4,p5)
    return err


def peval(x,p):
    A,B,C,D = p
    return logistic4(x,A,B,C,D)

def peval5(x,p):
    p1,p2,p3,p4,p5 = p
    return logistic5(x,p1,p2,p3,p4,p5)
'''
x = np.linspace(0,20,20)
p1,p2,p3,p4,p5 = 0.5,2.5,8,7.3,3.9
y_true = logistic5(x,p1,p2,p3,p4,p5)
y_meas = y_true + 0.2*npr.randn(len(x))

# Initial guess for parameters
p0 = [0,1,1,1]

# Fit equation using least squares optimization
#plsq = leastsq(residuals,p0,args=(y_meas,x))
plsq = least_squares(residuals5,p0,f_scale=0.1,args=(y_meas,x))
#plsq = least_squares(residuals5,p0,f_scale=0.1,args=(x,y_meas))

plt.plot(x,peval5(x,plsq.x),x,y_meas,'o',x,y_true)
plt.title('Least-squares 5PL fit noisy data')
plt.legend(['Fit','Noisy','True'],loc = 'upper left')
for i , (param, actual, est) in enumerate(zip('ABCDE',[p1,p2,p3,p4,p5],plsq.x)):
    plt.text(10,3-i*0.5,'%s = %.2f, est(%s) = %.2f' % (param, actual, param, est))
plt.savefig('logis_text5_new.png')
#plt.show()
'''
'''
x = np.linspace(0,20,20)
'''
x = [1,2,3,4,5,6,7,8,9,10,11,12,13]
x = np.array(x)
A,B,C,D = 0.5,2.5,8,7.3
y_true = logistic4(x,A,B,C,D)
y_meas = y_true + 0.2*npr.randn(len(x))
print type(x)
#x = [1,2,3,4,5,6,7,8,9,10,11,12,13]
# Initial guess for parameters
p0 = [0,1,1,1]

# Fit equation using least squares optimization
#plsq = leastsq(residuals,p0,args=(y_meas,x))
plsq = least_squares(residuals,p0,f_scale=0.1,args=(y_meas,x))


plt.plot(x,peval(x,plsq.x),x,y_meas,'o',x,y_true)
plt.title('Least-squares 4PL fit noisy data')
plt.legend(['Fit','Noisy','True'],loc = 'upper left')
for i , (param, actual, est) in enumerate(zip('ABCD',[A,B,C,D],plsq.x)):
    plt.text(10,3-i*0.5,'%s = %.2f, est(%s) = %.2f' % (param, actual, param, est))
plt.savefig('logis4_x_test.png')
#plt.show()


